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Abstract 

We propose extension of the algorithm for numerical modelling of Bose-Einstein correla- 
tions (BEC), which was presented some time ago in the literature. It is formulated on quantum 
statistical level for a single event and uses the fact that identical particles subjected to Bose 
statistics do bunch themselves, in a maximal possible way, in the same cells in phase-space. 
The bunching effect is in our case obtained in novel way allowing for broad applications and 
fast numerical calculations. First comparison with e"'"e~ annihilations data performed by us- 
ing simple cascade hadronization model is very encouraging. 
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Bose-Einstein correlations (BEC) between identical bosons are since long time of special interest 
because of their ability to provide space-time information about multiparticle production processes 
Q . This is particulary true in searches for a proper dynamical evolution of heavy ion collisions 
(QGP) 1^. However, such processes, because of their complexity, must be modelled by means of 
Monte Carlo event generators probabilistic structure of which prevents a priori the genuine 
BEC (which are of purely quantum statistical origin). One can only attempt to model BEC in 
some way aiming to reproduce the two-particle correlation function measured experimentally and 
defined, for example, as ratio of the two-particle distributions to the product of single-particle 
distributions: 

Ni{pi)Ni{pj) 

This is done either by changing the original output of these generators by artificially bunching in 
phase-space (in a suitable way) the finally produced identical particles [§, D or by constructing 
generator, which allows to account from the very beginning for the bosonic character of produced 
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secondaries [|6|[J. In former case the simplest approach is to shift (in each event) momenta of ad- 
jacent hke-charged particles in such a way as to get desired C2{Q) Q and to correct afterwards 
for the energy-momentum imbalance introduced this way. Much more physical is the method de- 
veloped in , which screens all events against the possible amount of bunching they are already 
showning and counts them as many times as necessary to get desirable C2(Q)^. The original 
energy-momentum balance remains in this case intact whereas the original single particle distri- 
butions are changed (this fact can be corrected by running again generator with suitably modified 
input parameters). The size parameters occuring in weights bear no direct resemblance to the size 
parameter R obtained by directly fitting data on C2{Q) in eq.(^ by, for example, simple gaussian 
form: C2{Q) = c [l + Acxp {—R^Q^)] (where c is normalization constant, A the so called chaoticity 
and R the size parameter). They represent instead the corresponding correlation lengths between 
the like particles [Q. 

The approach proposed in represents different philosophy of getting desired bunching. Here 
one uses specific generator, which groups (already on a single event level) bosonic particles of the 
same charge in a given cell in phase space according to Bose-Einstein distributiorrl 



PiE^) ~ exp 



Hi {p - Ei) 



(2) 



Here rii is their multiplicity and Ei energy, the energy-momentum and charge conservations are 
strictly imposed by means of the information theory concept of maximazing suitable information 
entropy. The parameters T and /i are therefore two lagrange multipliers with values fixed by the 
energy-momentum and charge conservation laws, respectively. Such distribution is typical example 
of nonstatistical fluctuations present in the hadronizing source. With only one additional param- 
eter 5y, which denotes the size of the elementary emitting cell in phase-space (in it means in 
rapidity), one gets at the same time both the correct BEC pattern (i.e., correlations) and fluctu- 
ations (as characterized by intermittency) . This is very strong advantage of this model, which 
is so far the only example of hadronization model, in which Bose-Einstein statistics is not only 
included from the very beginning on a single event level, but it is also properly used in getting final 
secondaries. In all other approaches at least one of the above elements is missing. The shortcom- 
ing of method are numerical difficulties to keep the energy-momentum conservation as exact as 
possible and limitation to the specific event generator only. 



In the present work we propose generalization of this approach making it applicable to other 
generators. Namely, following the same reasoning as in § , we propose different method of intro- 
ducing desired bunching. In Q the generator itself provided particles satisfying Bose statistics (in 
the sense mentioned above). In the general case one has to find the possible bosonic configurations 
of secondaries existing already among the produced particles. The point is that nonstatistical fluc- 
tuations present in each event generator result in a nonuniform (bunched) distributions of particles 
produced in a given event in momentum space. They resemble Bose distribution provided by the 
statistical event generator of the only difference being that particles in such bunches usually 
have different charges allocated to them by event generator, whereas in Q particles were of the 

^The specific approaches proposed for LUND model and the afterburner method discussed in which we 
shall not discussed here, belong to first cathegory. 

^Technically this is realised by multiplying each event by a special weight calculated using the output provided 
by event generator. 

^Similar concept of elementary emmiting cells has been also proposed in bl. 
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same charge. We propose therefore to maximally equalize charges of particles in such bunches (to 
the extent limited only by the overall charge conservation) superseeding the initial charge alloca- 
tion provided by event generator (keeping only intact the total number of particles of each charge 
it gives). This is supposed to be done in each single event. Both the original energy- momentum 
pattern provided by event generator and all inclusive single particle distributions are left intact^ 

It is instructive to look at this problem from yet another point of view. Namely, it can be 
perceived as an attempt (cf. [Q) to model correlations of fluctuations present in the system, as 
given by: 

(nin2) = (ni)(n2) + ((ni - (ni)) (n2 - (712))) 

= {ni){n2) + p(T{ni)(7(n2). (3) 

Here cr(n) is dispersion of the multiplicity distribution P{n) and p is the correlation coefficient 
depending on the type of particles produced: p = +1,-1,0 for bosons, fermions and Boltzmann 
statistics, respectively. The proposed algorithm should provide us with C2{Q), which is in fact a 
measure of correlation of fluctuations because 

CiQ - \P.-P,\) ^ ^-f^l^^^f ^ 1 + ^7^7^- (4) 

(rii [Pi ) ) {rij [pj ) ) {n, [p, ) ) {uj (pj ) ) 

To get p > (in which we are interested here) it is enough to select one of the produced particles, 
allocate to it some charge, and then allocate (in some prescribe way) the same charge to as many 
particles located near it in the phase space as possible (limited only by the charge conservation 
constraint^). In this way one forms a cell in phase-space, which is occupied by particles of the same 
charge only. This process should then be repeated until all particles are used and it should be done 
in such way as to get geometrical (Bose-Einstein) distribution of particles in a given cell. We stress 
again that this procedure does not alter neither the original energy-momentum distributions nor 
the spatio-temporal pattern of particles provided by event generator. It only changes the charge 
flow pattern it provides retaining, however, both the initial charge of the system and its total 
multiplicity distribution. Therefore this method works only when we can resign from controlling 
the charge flow during hadronization process. 



The procedure of formation of such cells is controlled by a parameter P being the probability 
that given neighbor of the initially selected particle should be counted as another member of the 
newly created emmiting cell in phase space. Notice that any selection procedure leading to a 
geometrical particle distribution in cells (in which case cr n >) results in maximization of the 
second term in the eq. (^ . In particular it can be realized by the following algorithm of allocation 
of charges (the ni = nl^' -I- n[ -I- n[''^ is the number of particles of different charges provided by 
our event generator in the l*^ event, {pj} and {xi} are, respectively, their energy-momenta and 
spatio-temporal positions, which we keep intact): 

^The necessary changes to event generator used introduced by such procedure will be discussed later. What 
we propose here is to resign from the part of the information provided by event generator concerning the charge 
allocation to produced particles. This can be regarded as introduction of quantum mechanical element of uncertainty 
to the otherwise classical scheme of generator used (however, it differs completely from the usual attempts to 
introduce quantum mechanical effects discussed in |lo|). 

^The energy-momentum constraint is taken caxe by the generator itself and is not affected by our algorithm. 
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(1) The SIGNis chosen randomly from: "+", "-"or "0" pool, with probabilities given by p\ 
n^^^ /ni, p\ = n\ ^ jrii and = nf^ jni. It is attached to particle (i) chosen randomly 
from the particles produced in this event and not yet reassigned new charges. 

(2) Distances in momenta, Sij{p) = \pi — Pj\, between the chosen particle (i) and all other par- 
ticles (j) still without signs are calculated and arranged in ascending order with j = 1 
denoting the nearest neighbor of particle (i). To each Sij{p) one assignes some probability 

P{i,j) e (0,1). 

(3) A random number r G (0, 1) is selected from a uniform distribution. If nf^'-'^ > 0, i.e., if 
there are still particles of given SIGN with not reassigned charges, one checks the particles (j) 
in ascending order of j and if r < P{i, j) then charge SIGNis assigned also to the particle (j), 
the original multiplicity of particles with this SIGNis reduced by one, nf^'-^^ = nf^'^^ — 1, 
and the next particle is selected: (j) {j + 1). However, if r > P{i,j) or nf^^^ — then 
one returns to point (1) with the updated values of p[ ' and pi"'- Procedure finishes 
when nl^"* — n| •* = = 0, in which case one proceeds to the next event. 

As can be easily checked, this algorithm results in geometrical (Bosc-Einstein) distribution of par- 
ticles in the phase-space cells formed by our procedure (with mean multiplicity P/{1 — P) for 
constant P(ij) = P case) accounting therefore for their bosonic character (i.e., for Bose-Einstein 
statistics they should obey)|. 

We shall illustrate now action of our algorithm on simple cascade model of hadronization 
(CAS) (in its one-dimensional versions and assuming, for simplicity, that only direct pions are 
produced) . In CAS the initial mass M hadronizes by series of well defined (albeit random) 
branchings [M — > Mi -I- M2) and is endowed with a simple spatio-temporal pattern. It shows no 
traces of Bose-Einstein statistics whatsoever. However, as can be seen in Fig. 1, when endowed 
with charge selection provided by our algorithm, a clear BEC pattern emerges in C2 (Q) . Two kind 
of choices of probabilities are shown in Fig. 1. First is constant P — 0.75 and P = 0.5. It leads 
to a pure geometrical distribution of number of particles allocated to a given cell and corresponds 
to situation already encountered in Its actual value is so far a free parameter replacing, in a 
sense, the parameter Sy in However, whereas in the size of emitting cells was fixed, in our 
case it is fluctuating. The other is what we call the "minimal" weight constructed from the output 
information provided by CAS event generator: 



P{ij) = exp 



(5) 



where Sij (x) = \xi — Xj \ and Sij (p) — \pi —Pj\- In this way one connects P with details of hadroniza- 
tion process by introducing to it a kind of overlap between particles as a measure of probability of 
their bunching in a given emitting cell. 



As we have checked out the BEC effect occuring here depends only on the (mean) number of 
particles of the same charge in phase-space cell and on the (mean) numbers of such cells. This 

^It is important to realize that, because we do not restrict a priori the number of particles which can be put in a 
given cell, we are automatically getting BEC of all orders (even if we use only two particle checking procedure at a 
given step in our algorithm). It means that C2(Q = 0) calculated in such environment of the possible multiparticle 
BEC can exceed 2 (cf. 
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depends on P, the bigger P the more particles and bigger C2{Q — 0); smaller P leads to the 
increasing number of cells, which, in turn, results in decreasing C2{Q = 0), as already noticed in 
1^ . For small energies the number of cells decreases in natural way while their occupation remains 
the same (because P is the same), therefore the corresponding 6*2(0) is bigger, as seen in Fig. 1. 
The fact that there is tendency to have 6*2(0) > 2 for larger P means that one has in this case 
more cells with more than 2 particles allocated to them, i.e., it is caused by the influence of higher 
order BEC. Therefore the "sizes" R obtained from the exponential fits to results in Fig. 1 (like 
C2{Q) ~ 1 + A • exp(— Q • R) where A being usually called chaoticity parameter Q) correspond to 
the sizes of the respective elementary cells rather than to sizes of the whole hadronizing sources 
itself. For P = 0.5 the "size" R varies weakly between 0.66 to 0.87 fm from M = 10 to 100 GeV 
whereas for the "minimal" weight (^ it varies from 0.64 to 0.44 fniQ 

So far we were considering only single sources. Suppose now that source of mass M consists 
of a number (n; — 2^) of subsources hadronizing independently. It turns out that the resulting 
C2's are very sensitive to whether in this case one applies our algorithm of assigning charges to 
all particles from subsources taken together ("Split" type of sources) or to each of the subsource 
independently ("Indep" type of sources), cf. Fig. 2. Whereas the later case results in the similar 
"sizes" R (defined as before) with C2(Q = 0) — 1 = A falling dramatically with increasing k (roughly 
like 1/2*^, i.e., inversely with the number of subsources, ni, as expected fom 1^), the former leads 
to roughly the same C2{Q ~ 0) but the "size" R is now increasing substantially being equal to, 
respectively, 0.87 fm, 1.29 fm, 1.99 fm and 3.35 fm for P = 0.5 and 0.57 fm, 3.26 fm, 4.01 and 5.59 
fm for the "minimal" weight (|^). Fig. 3 contains example of our "best fit" to the e+e" annihilation 
data on BEC by DELPHI Collaboration ||l|l for M = 91.3 GeV (which can be obtained only for 
two or three subsources, as shown there). 

To summarize: we propose a new and simple method of numerical modelling of BEC. It is based 
on reassigning charges of produced particles in such a way as to make them look like particles satis- 
fying Bose statistics, conserves the energy-momenta and does not alter the spatio-temporal pattern 
of events or any single particle inclusive distribution (but it can change the distributions of, sepa- 
rately, charged and neutral particles leaving, however, the total distribution intact). It is intended 
to be a kind of suitable extension of the algorithm presented in |^ , such that can be applied to 
essential any event generator in which such reassignment of charges is possible. It amounts to the 
changes in physical picture of the original generator. The example of CAS is very illustrative in 
this respect. In it, at each branching vertex one has, in addition to the energy-momentum conser- 
vation, imposed strict charge conservation and one assumes that only (0) (H — ), (+) (+0) 
and (— ) (—0) transitions are possible. It means that there are no multicharged vertices (i.e., 
vertices with multiple charges of the same sign) in the model. However, after applying to the 
finally produced particles our charge reassignment algorithm one finds, when working the branch- 
ing tree "backwards", that precisely such vertices occur now (with charges "(-I--I-)", or "(- -)", for 
example). The total charge is, however, still conserved as are the charges in decaying vertices (i.e., 
no spurious charge is being produced because of our algorithm). It is plausible therefore that to 
get BEC in an event generators it is enough to allow for cumulation of charges of the same sign 
at some points of hadronization procedure modelled by this generator. This would lead, however, 

'^This should be contrasted with the "real" (mean) sizes of CAS sources changing from 0.43 fm for M = 10 GeV 
to 1.82 fm for M = 100 GeV Q. 
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to extremely difficult numerical problem with ending such algorithms without producing spurious 
multicharged particles not observed in nature^. 

We find that value of C2{Q = 0) (defining chaoticity parameter A) depends inversely on the 
number of elementary cells, in the way already discussed in and that "radius" R extracted 
from the exponential fits is practically independent on the size of the source, provided it is a single 
one. In the case when it is composed of a number of elementary sources, R increases with their 
number, unless they are treated independently by our algorithm. It is because one has in this case 
a higher density of particles. This results in smaller average Q, and this in turn leads to bigger 
On the contrary, for the independently treated sources the density of particles subjected to 
our algorithm does not change, hence the average Q and R remain essentially the same. How- 
ever, because in this case the infiuence of pairs of particles from different subsources increases, the 
effective A = 6*2(0) — 1 now decreases substantially (as was already observed in |^). It should 
be mentioned at this point that our "Indep" type sources can probably be used as a possible ex- 
planation of the so called inter- BEG problem, i.e., the fact that essentially no BEG is being 
observed between pions originating from a different W in fully W~^W~ final states [0. It can be 
understood by assuming that produced W^'s should be treated as "Indep" type sources for which A 
falls dramatically. Finally, we would like to stress that our algorithm leads to strong intermittency 
showing up after its application. It means that with such algorithm (which is very efficient and 
fast for all multiplicities) we can already attempt to fit experimental data by applying it to some 
more sophisticated fragmentation schemes than that provided by GAS. This will be done elsewhere. 
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Figure 1: Examples of BEC patterns obtained for M = 10, 40 and 100 GeV for constant weights 
P = 0.75 (stars) and P = 0.5 (full symbols) and for the weight given by eq.(||) (open symbols). 
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Figure 2: Examples of BEC for different number of subsources (n; = 2*'', k — 1,2,3 existing in 
the source M — 100 GeV. Left panel is for the "Split" and right panel for the "Indep" types of 
sources, as discussed in text. 
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Figure 3: The examples of the "best fits" to the e+e annihilation data on BEC by DELPHI ||lJ 
using 2 (a) or 3 (6) subsources with P — 0.23 and P = 0.3, respectively. 
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